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DYNAMICS  OF  AN  UNSTEADY  DIFFUSION  FLAME:  EFFECTS 
OF  HEAT  RELEASE  AND  GRAVITY 


INTRODUCTION 

Experiments  on  laminar  diffusion  flames  have  shown  that  gravity  affects  the  flame 
length  and  width  as  well  as  its  extinction  characteristics  (1-4).  These  studies  have 
been  conducted  in  drop  towers  and  have  focused  on  fuel  jets  with  very  low  velocities 
of  less  than  50  cm/s.  Although  these  experiments  have  increased  our  basic  under¬ 
standing  of  laminar  diffusion  flames  by  emphasizing  the  importance  of  bouyancy,  it 
is  not  clear  how  to  apply  these  results  to  higher- velocity  flames  which  are  unsteady 
or  fluctuating.  Studying  higher- velocity  fuel  jets  from  larger  nozzles  is  more  difficult 
experimentally  because  the  flames  can  be  quite  long  and  the  instabilities  may  not 
have  time  to  evolve  during  a  single  experiment.  Through  numerical  simulations,  we 
can  examine  an  unsteady  flame  with  and  without  gravity  in  the  kind  of  detail  that 
is  not  practical  in  an  experiment. 

Two  types  of  instabilities  axe  observed  in  low-speed  diffusion  flames  (5,6).  The 
high-frequency  structures  grow  from  Kelvin- Helmholtz  instabilities  at  the  interface 
between  the  high-velocity  and  low-velocity  fluid  and  typically  have  frequencies  of 
a  few  hundred  Hertz.  The  low-frequency  structures  form  in  the  region  outside  the 
flame  zone  with  typical  frequencies  of  10-20  Hertz. 

This  paper  examines  the  effect  of  heat  release  and  gravity  on  the  formation  and 
evolution  of  these  two  types  of  instabilities  by  presenting  a  series  of  time-dependent, 
two-dimensional  simulations  of  an  axisymmetric  H2-N2  jet  in  a  coflowing  air  stream. 
The  calculations  include  convection,  thermal  conduction,  molecular  diffusion,  vis¬ 
cosity,  chemical  reactions  with  energy  release,  and  gravitational  forces.  This  model 
is  based  on  the  one  developed  by  Laskey  (7),  which  includes  a  new  algorithm  for 
convective  transport  developed  by  Patnaik  et  al.  (8).  Previously,  Laskey  (9)  pre¬ 
sented  computations  of  diffusion  flames  of  the  type  presented  here  and  Patnaik  et 
al.  (10)  used  a  similar  model  to  study  the  stability  properties  of  very  low-speed 
premixed  flames.  All  of  these  efforts  have  tested  the  various  parts  of  the  model  and 
have  given  credibility  to  its  overall  validity. 

Manuscript  approved  May  25,  1990. 
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These  computations  are  different  from  previous  ones  reported  by  Laskey  et 
al.  (9)  and  Ellzey  et  al.  (11)  for  several  reasons.  Greater  resolution  and  improved 
boundary  conditions  now  allow  correct  zero-gravity  computations.  The  energy- 
release  model  now  properly  limits  the  final  temperatures  allowed  and  no  longer 
produces  a  strong  recirculation  zone  at  the  jet  exit.  Viscosity  has  been  shown  to  be 
very  important  in  diffusion  flames  at  these  velocities  (11)  and  is,  therefore,  included 
in  all  of  the  calculations  presented  in  this  paper. 

NUMERICAL  METHOD 

The  numerical  model  consists  of  separate  algorithms  for  the  various  processes,  and 
these  algorithms  are  coupled  by  timestep  splitting  methods.  Table  I  is  an  outline 
of  one  computatonal  timestep.  Given  a  set  of  initial  values  for  the  basic  vari¬ 
ables,  an  approriate  computational  timestep  is  estimated  based  on  accuracy  and 
stability  criteria.  Then  the  effects  of  thermal  conduction  are  evaluated  using  a  two- 
dimensional  explicit  finite-difference  model  (7).  Thermal  conductivities,  for  the 
individual  species  were  calculated  from  kinetic  theory  over  the  temperature  range 
300  to  3300  K,  these  values  were  fit  to  a  third-order  polynomial,  and  then  are  used 
to  calculate  the  mixture  thermal  conductivity  (13).  Molecular  diffusion  is  included 
using  an  explicit  finite-difference  formulation.  First,  the  diffusion  velocities  are  cal¬ 
culated  according  to  Fick’s  law  and  then  corrected  (13)  to  satisfy  the  requirement 
that  the  sum  of  the  diffusion  fluxes  is  zero.  Binary  diffusion  coefficients,  calculated 
from  kinetic  theory  (14),  are  used  to  compute  the  diffusion  coefficients  for  a  partic¬ 
ular  species  in  a  mixture  (13).  The  viscosity  coefficients  fit,  calculated  from  kinetic 
theory  over  the  temperature  range  300  to  3000K  and  fit  to  a  third  order  polynomial, 
were  used  to  compute  the  mixture  viscosity  (15).  The  model  for  chemical  reactions 
and  heat  release  is  an  extension  of  the  Parametric  Diffusion  Reaction  Model  (7,12), 
which  is  designed  to  replace  the  integration  of  the  full,  detailed  set  of  ordinary 
differential  equations  representing  the  chemical  kinetics.  A  single,  global  reaction 
is  used  but  the  reaction  is  not  instantaneous.  Instead,  the  finite  reaction  rate  is 
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Table  I.  One  Timestep  in  the  Diffusion-Flame  Model 


Given  Initial  Variables 

1.  Determine  Timestep 

2.  Thermal  Conduction 

Integrate  from  t  to  t  4-  At: 

Calculate  Aei.  Do  not  update  any  variables.  Subcycle  as  necessary. 

3.  Ordinary  Diffusion 

Integrate  from  t  to  t  +  At: 

Only  update  {n,(z)}.  Calculate  Ae2-  Subcycle  as  necessary. 

4.  Viscosity 

Integrate  from  t  to  t  +  At: 

Only  update  pv.  Calculate  Aea. 

5.  Chemical  Reactions 

Integrate  from  t  to  t  -f  At: 

Only  update  {rij(x)}.  Calculate  Ae.4. 

6.  Convective  Transport 

Integrate  from  t  to  t  +  At: 
x  direction  transport,  then  update  p,  pv,E,rii. 
y  direction  transport,  then  u  update  p,  pv,E,rii. 

Implicit  correction,  then  update  p ,  e,  and  E. 

7.  Increment  Time  and  go  to  1. 


determined  such  that  the  maximum  temperature  in  a  one-dimensional  transient 
diffusion  flame  is  the  adiabatic  flame  temperature  for  a  stoichiometric  mixture  of 
the  fuel  and  oxidizer.  The  transport  of  density,  momentum,  energy,  and  individual 
species  density  is  accomplished  through  the  high-order  implicit  method,  BIC-FCT 
(8).  This  involves  an  explicit  step,  based  on  the  standard  FCT  algorithm  (16),  and 
then  an  implicit  correction. 

The  general  timestep  splitting  approach  for  coupling  the  various  physical  pro¬ 
cesses  was  developed  for  slow-flow  implicit  calculations.  In  these  computations,  the 
change  in  internal  energy  resulting  from  each  individual  process  is  not  incorporated 
into  the  solution  as  soon  as  it  is  computed,  but  instead  is  accumulated,  as  indicated 
by  the  {e*}  in  Table  1.  The  entire  change  in  internal  energy  is  then  added  to  the 
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energy  equation  in  the  fluid  convection  step  6.  The  coupling  technique  has  been 
described  by  Oran  and  Boris  (16),  and  a  modification  by  Patnaik  et  al.  (17)  has 
been  shown  to  allow  for  a  greater  addition  of  energy  per  timestep  while  maintaining 
numerical  stability. 

In  essence,  the  model  solves  the  time- dependent  two-dimensional  conservation 
equations  for  mass  density,  p,  momentum,  pv,  and  total  energy,  E  and  these  are 
coupled  to  models  for  chemical  reactions  among  the  species  {n,}  with  subsequent 
heat  release,  molecular  diffusion,  thermal  conduction,  viscosity,  and  gravitatonal 
forces.  Additional  equations  include  the  perfect  gas  equation  of  state  and  a  relation 
between  the  internal  energy  and  the  pressure.  The  specific  set  of  equations  and 
more  detailed  discussions  of  the  numerical  methods  are  given  in  References  (7)  and 
(12). 

The  computations  described  in  this  paper,  using  the  enlarged  computational 
grid  and  including  all  of  the  physical  processes,  require  0.7  s/computational  timestep 
on  a  Cray  YMP.  This  means  that  a  typical  calculation,  about  50,000  timesteps, 
requires  about  10  hours  of  computer  time. 

APPLICATION  TO  UNSTEADY  DIFFUSION  FLAMES 

The  computational  grid  for  the  region  near  the  jet  and  the  initial  conditions  are 
shown  in  Figure  1.  The  full  domain  is  10  cm  x  172  cm  and  consists  of  128  x  224  cells. 
Cells  of  approximately  0.02  cm  are  concentrated  around  the  jet  exit.  Beginning  at 
r  =  1  cm,  the  size  of  each  cell  is  increased  by  3%  over  the  size  of  its  neighboring  cell 
for  all  simulations.  The  cells  in  the  axial  direction  for  all  simulations  are  stretched 
by  3%  starting  at  z  =  1  cm.  A  fuel  mixture  consisting  of  78%  H2  and  22%  AT2 
by  volume  flows  through  a  jet  of  radius  0.5  cm  at  10  m/s  at  the  lower  boundary. 
Air  flows  through  the  outer  annular  region  between  r  =  0.5  and  r  =  10.0  cm  at  30 
cm/s.  The  outer  boundary  at  r  =  10.0  cm  is  a  free-slip  wall.  The  inner  boundary 
at  r  =  0.0  is  the  jet  centerline.  An  outflow  boundary  is  specified  at  z  =  172  cm 
where  the  pressure  is  adjusted  to  atmospheric. 
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RESULTS 


Nonreacting  Jet 

Figure  2  shows  the  instantaneous  contours  of  axial  and  radial  velocity  and  mole 
fraction  late  in  the  simulation  of  zero-gravity  nonreacting  jet.  Kelvin-Helmholtz 
instabilities  occur  near  the  jet  exit  leading  to  vortical  structures  that  then  con- 
vect  downstream.  These  structures,  which  transport  fuel  and  oxidizer  radially  and 
broaden  the  mixing  zone,  weaken  substantially  in  the  first  ten  jet  diameters.  Small 
radial  velocities,  not  evident  in  the  contours,  still  exist  at  this  point.  Figure  3 
shows  the  mean  and  rms  velocity  for  the  nonreacting  jet  at  three  axial  locations. 
At  z  =  0.5  cm,  the  mixing  region  is  narrow  with  only  small  fluctuations  of  a  few 
cm/s.  At  z  =  1.0  cm,  the  instabilities  result  in  large  fluctuations  across  the  entire 
jet  core.  By  z  =  10  cm,  there  are  small  fluctuations  across  the  entire  jet  region. 
The  results  for  the  nonreacting  jet  with  gravity  are  not  distinguishable  from  those 
for  the  same  jet  in  zero  gravity,  and  so  axe  not  shown  here. 

Reacting  Jet 

Instantaneous  contours  late  in  the  calculation  of  the  reacting  jet  in  zero  gravity, 
Figure  4,  show  that  the  volumetric  expansion  and  the  change  in  temperature  have 
a  significant  effect  on  the  flow.  The  radial  velocities  arise  from  the  expansion  at  the 
flame  front  but  are  relatively  uniform.  The  axial  velocity  and  concentration  fields 
are  steady  in  time.  Figure  5,  the  mean  and  rms  velocity  for  this  case,  show  that 
the  mixing  region  is  wider  due  to  the  expansion.  Fluctuations  are  insignificant  and 
not  visible  on  the  plot. 

Figure  6  shows  that  gravity  changes  the  flow  significantly.  Instabilities  form 
outside  the  reaction  zone  in  the  region  with  large  temperature  and  density  gradients. 
The  maximum  radial  velocity  is  approximately  30  cm/s  and  occurs  at  the  center 
of  the  structure.  The  concentration  and  temperature  fields  are  distorted  as  these 
instabilities  convect  downstream.  The  flame  front  lies  at  the  fuel-oxidizer  interface 
in  the  region  of  maximum  temperature  and  fluctuates  in  time.  Figure  7  shows  a 


5 


time  sequence  of  the  H2O  mole- fraction  contours.  In  the  first  frame,  a  bulge  is 
developing  on  the  outside  of  the  H2O  contours.  In  subsequent  frames,  it  rolls  up 
and  moves  downstream.  In  the  final  frame,  it  is  moving  out  of  the  domain  shown  as 
the  next  instability  forms  below  it.  These  outer,  slower-moving  vortical  structures 
occur  at  approximately  15  Hz. 

DISCUSSION  AND  CONCLUSIONS 

Comparisons  of  the  four  computations  of  the  10  m/s  H2  —  N2  jet  into  the  SO  cm/s 
coflowing  air  background  shows  that  gravity  and  heat  release  interact  substantially 
to  change  the  flow.  Without  chemical  reactions  and  subsequent  heat  release,  gravity 
does  not  noticeably  change  the  velocity  or  concentation  fields.  Even  though  there 
are  significant  density  gradients  between  the  H2  —  N2  fuel  jet  and  the  co-flowing 
air,  these  gradients  occur  in  a  region  of  relatively  high  velocity  where  momentum 
effects  dominate. 

In  the  reacting  jet,  there  are  significant  density  gradients  in  the  coflow  region 
where  the  velocity  is  low.  These  gradients  are  due  to  the  conduction  of  heat  away 
from  the  reaction  zone.  In  this  region,  the  bouyant  forces  dominate  and  large 
instabilities  form.  These  have  been  observed  in  experiments  for  many  years  (6, 
18-20)  and  are  considered  to  be  responsible  for  flame  flicker. 

Volumetric  expansion  and  the  effects  of  changing  temperature  stabilize  the  mix¬ 
ing  region  of  the  reacting  jet.  The  increase  in  viscosity  with  temperature  accounts 
for  part  of  the  stabilization  but  analytical  results  show  that  inviscid  instabilities  are 
also  damped  by  heat  release  (21).  Preliminary  computations  with  constant  viscosity 
indicate  that  the  stabilization  effect  due  to  the  change  in  viscosity  with  temperature 
may  be  insignificant  compared  to  the  effect  of  heat  release.  Previous  calculations 
(11)  show  that  even  without  including  viscosity,  heat  release  reduces  the  strength 
of  the  Kelvin-Helmholtz  instability. 

Future  computations  are  proceeding  in  several  different  directions.  First,  we 
are  considering  the  downward-propagating  diffusion  flame  and  how  this  differs  from 
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upward  and  zero-gravity  flames.  Second,  we  are  x«.  iucing  the  coflow  velocity  so  that 
the  computations  have  the  same  parameters  as  recent  experiments  at  the  Air  Force 
Wright  Aeronautical  Laboratory.  At  that  point,  detailed  comparisons  will  be  made 
between  the  computations  and  experimental  results.  We  are  continuing  to  develop 
the  energy-release  model  so  that  the  energy  release  as  a  function  of  temperature 
is  better  represented.  Finally,  we  have  been  investigating  new  types  of  computers 
that  might  allow  full-chemistry  or  three-dimensional  computations  of  such  flames. 
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Figure  1.  Computational  domain  and  initial  conditions  for  the  compu¬ 
tations  of  a  Hj  -  Ni  jet  into  coflowing  air.  Note  that  the  figures  only 
show  the  part  of  the  full  domain  with  the  high  resolution. 


Figure  2.  Contours  of  (a)  radial  velocity,  (b)  axial  velocity,  (c)  mole  fraction  Hi,  (d)  mole 
fraction  02  for  a  nonreacting,  zero-gravity  jet  of  H2  -  Ni  into  coflowing  air.  Dimensions 
are  in  cm,  velocities  are  in  cm/s. 
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Figure  3.  Mean  and  rms  velocity  for  the  zero-gravity  nonreacting  jet  at  three  axial  locations. 
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Figure  4.  Contours  of  (a)  radial  velocity,  (b)  axial  velocity,  (c)  mole  fraction  Hit  (d  >  mole  fraction  ()  .,  (d)  temperature 
for  a  zero-gravity  diffusion  flame  formed  between  a  H2  -  N2  jet  and  coflowing  air.  Dimensions  are  in  cm,  velocities 
are  in  cm/s,  temperature  is  in  K. 
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Figure  5.  Mean  and  rms  velocity  for  the  zero-gravity  diffusion  dame  at  three  axial  locations. 
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Figure  7.  Sequence  of  contours  of  H?0  mole  fraction  showing  the  formation  of  the  large  vortical  structures  in  the 
coflow  region.  The  time  interval  between  each  frame  iB  11.25  ms. 


